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Abstract 

A generalised integer S Ising spin glass model is analysed using the replica 
formalism. The bilinear couplings are assumed to have a Gaussian distribu- 
tion with ferromagnetic mean < Jjj >= Jq. Incorporation of a quadrupolar 
interaction term and a chemical potential leads to a richer phase diagram 
with transitions of first and second order. The first order transition may be 
interpreted as a phase separation, and contrary to what has been argued previ- 
ously, it persists in the presence of disorder. Finally, the stability of the replica 
symmetric solution with respect to fluctuations in replica space is analysed, 
and the transition lines are obtained both analytically and numerically. 
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I. INTRODUCTION 



The much studied SK model of spin glasses may be generalised in different ways. The 
standard spin glass model, given by 

H = - ^2 Jij s i s 3 witn s i = ±1) (1) 

<ij> 

may be extended for instance by allowing values s = 0, ±1, ±2, . . . , ±S for the spin variables. 
It is then possible to consider higher order interactions such as KJ^yS^s 2 , or a chemical 
potential such as A J2i s 1- Such generalisations can be regarded as extensions of the Blume- 
Emery-Griffiths model (BEG) [|J. Indeed the BEG model allows s = 0, ±1 and takes into 
account the aforementioned higher order interaction and chemical potential. However, in 
the standard formulation, the bilinear couplings are neither frustrated nor disordered: they 
are ferromagnetic. The s = degree of freedom has been used to model a diluted lattice gas 
and leads to a first order transition separating magnetic and non-magnetic phases. Thus 
our generalisation of the BEG Hamiltonian is a good way to study the influence of disorder 
and frustration on first order transitions. 

The BEG model has been studied in various contexts and a mean field approximation 
has been given by Blume, Emery and Griffiths 0. They introduced the model in order 
to study multicritical phenomena associated with physical systems such as binary mixtures. 
For an overview, see the review article of Lawrie and Sarbach ||18|| . To improve the mean field 



results and to cope with the finite dimensionality of real physical systems, several different 
techniques have been applied to the BEG model, e.g., Tucker |24| applied the exponential 
operator technique of Honmura and Kaneyoshi fllfj to the isotropic BEG, and Fittipaldi et 



al. [12fl applied it to the anisotropic BEG model. In addition the model has been treated 



in the cluster approximation by Tucker [25] and in the local mean field approximation by 
Maritan et al. 01. 



Several disordered BEG models have been studied. Berker et al. [Ill looked at a bond 



disordered BEG model and Ez-Zahraouy [M looked at a bond diluted BEG model. To our 
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knowledge no studies of the case of quenched disorder with frustration have been carried 
out for the BEG model, though Arenzon et al. [|TJ have considered a frustrated lattice gas 
model similar to the BEG model. The closest model to a generalised BEG studied in the 
literature was first analysed by Ghatak and Sherrington ||13|| . In their model, they considered 
s = 0, ±1, and the influence of a chemical potential, but they had no quadrupolar interaction 
and the ferromagnetic mean Jo of their bilinear coupling was assumed to vanish. The 
interesting point we wish to stress here is that the first order transition of the BEG model, 
which can be interpreted as a phase separation transition, persists in the GS generalisation. 
This is in contrast to the findings of Berker et al. which will be discussed below. 

In the present paper we generalise the standard SK spin glass Hamiltonian to incorporate 
the chemical potential, the quadrupolar interaction, randomly distributed bilinear couplings 
with a non-zero ferromagnetic mean, and other integer values for the spin variables. We 
restrict our study to the s = 0, ±1 case, which should cover in a qualitative way the general 



integer S'-spin models [ FT] . Furthermore we do not study the effect of disorder in the 
quadrupolar couplings; such effects have been considered in other models by Snowman et al. 
|20|| . The paper is organised as follows. In section 2 we introduce the model. In section 3 the 
free energy is derived in the replica symmetric approximation and in section 4 we give the 
model's phase diagram. In section 5 we discuss our results in the light of recent arguments 
of Berker et al. concerning the influence of quenched disorder on phase transitions of first 
order. In section 6 we analyse the stability of the replica symmetric solution and derive the 
lines of instability. Conclusions are drawn in the final section. 

II. THE MODEL 

We consider the model described by the BEG-SK Hamiltonian, where the spin variables 
are allowed to assume the values Sj = 0, ±1 

H = -\H Jaw + A £ *? - \k E ^ (2) 

<ij> i <ij> i 
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This model will be treated in the mean field approximation, i.e., in the infinite range 
limit. For the sake of simplicity we assume, as in the SK model, bilinear couplings with a 
Gaussian distribution about a non-zero mean < J^ >= Jo > allowing for ferromagnetic 
ordering. To avoid additional complexity of the model we consider only the case of positive 
quadrupolar coupling, K > 0, and a ferromagnetic mean J of the bilinear couplings. The 
chemical potential A is allowed to assume positive as well as negative values. 

The couplings must be rescaled for the present case of infinite ranged couplings in order 
to give a sensible free energy in the thermodynamic limit: J — > J/y/N, Jo — > Jq/N and 
K — > K/N, respectively. 

We consider now some special cases, which have already been treated in the literature. 
When J = and Jo > 0, we recover the non- frustrated and non-disordered standard BEG 
model, which describes a binary mixture (e.g., He 3 and He 4 ). When J = and K = 0, we 
recover the Gathak and Sherrington model. In order to cope with the problem of averaging 
over quenched disorder, it is necessary to average the free energy over the bond distribution 
Jif f = F {Jij} /N j ... We accomplish this by using the replica technique of Edwards 
and Anderson |J. This technique relies on the identity ln[Z] = lim n _»o -{Z n — 1); Z n is 
interpreted as the partition function of a n-fold replicated system s$ — > sf, a = 1, . . . , n. 
The average free energy may be computed using the prescription: 

/3/ = -lim -(Z^-l). 

We follow this standard procedure to average the logarithm of the partition function. 
The average of the n-fold replicated partition function over the disorder gives 

I ij - a ij 

+^ E E ^ 5>«) 2 +(3hJ2J2sA (3) 
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The Hubbard-Stratonowich transformation e 2 ^^' 5 "*' 2 = JVt e^^^" with the con- 

-1*2 

vention T>t = dt gives for the free energy, within the framework of the replica method 
at the saddle point, 

Pfn = b-h Y.< + \P 2 J 2 E 4 + \PK £ zl + V J 2 £ g 2 , - In [Z^ , (4) 
where the effective Hamiltonian and its partition function is given by 

#cfr = E1-/3A + (3Kz a + \p 2 J 2 z a }sl + Y}P J v m <* + P h K + £ g^s^ 

Zcs = tr Sa=0 ,±i,...,±sexp{i? eff [s Q ]} 

The quantities introduced by these transformations acquire the meaning of order param- 
eters: 

m a = < s a > 

1 - x a = z a = < s 2 a > (5) 

where the average is with respect to the effective Hamiltonian. These results are in fact 
valid for general integer spin values, but in what follows we restrict ourselves to the case 
of S = 1. In order to solve this model it is necessary to make assumptions on the order 
parameter matrix q a p and to propose an Ansatz. 

III. FREE ENERGY OF THE MODEL IN RS 

We limit ourselves to the simplest Ansatz for the order parameter matrix, i.e., we assume 
symmetry with respect to permutations of any pair of the replicas: q a p = q, Va ^ (3. The 
single indexed quantities are assumed to be independent of the replica index: m a = m and 
z a = z, Va. This leads to: 

Pf = \pJom 2 + \pKz 2 + 1/3 VV - JV (6) 
- J Vy In [l + 2 cosh [(3 J m + (3h + (3J^qy\ e ^z-pA+^j^-^j^ q 
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This equation reduces to the one analysed by Ghatak and Sherrington |j!3| , Lage and de 
Almeida |T7| and by Mottishaw and Sherrington ||T9[ if J = and if there is no quadrupolar 
coupling. In order to simplify the notation hereafter, we define 

Mv) = Tr^rtr s =o, ± i,.., ± ^e-^>]. (7) 

The effective Hamiltonian and its partition function are the replica symmetric equivalents 
of those defined earlier and obtained by one further Hubbard-Stratonowich transformation. 
The set of functions 4>k(y) reduces for the S = 1 model to just two functions: 

2cosh[f3(J m + yJ y /q + h)\ 

My) = ~7T* — tt^-To — 2 io r k even, 

ePA+tfPJtq-PKz + 2 cos h[/3( J m + yj^q + h)} 

M y)= ^mJom + yjy-^h)] forfcodd 

^A+ipjtq-pKz + 2 cos h[/3( J m + y J^g + /i)] 

Introducing an effective temperature-dependent quadrupolar coupling k = K + \f3J 2 , 

the mean field saddle point equations in the replica symmetric approximation are: 

df r 

— = J m - J T>y = 0, 

urn J 

^ = kz - k Jv y My) = o, (8) 

^ = ~PJ 2 q + \PJ 2 J Vy [My)} 2 = 0. 



IV. PHASE DIAGRAM 

Mean field phase diagrams are obtained by extremising the free energy with respect to 
the order parameters. They have been determined for the BEG-SK model by extremising 
numerically the free energy with respect to the order parameters z, m and q for any 
temperature T and chemical potential A. Following a line of constant chemical potential 
while varying the temperature — as is shown in figure [5] for two different but fixed values 
of the chemical potential — will occasionally reveal the onset of ordering. The examples 
include the appearence of non- vanishing values of the order parameters m or q and discon- 
tinous changes of the order parameter z which indicate a phase separation. We infer the 
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transition lines by monitoring the magnetisation m, the spin glass order parameter q, and 
the concentration z for fixed chemical potential A and varying the temperature T. This 
was performed for a range of different chemical potentials, sufficiently large to exhibit the 
different phenomenae. 

To present the numerical results we will make use of two commonly used phase diagram 
sections, the section of chemical potential versus temperature (A-T) and the section of 
concentration versus temperature (x-T). For the sake of completeness we first give the 
results for the non-disordered BEG model. Since in the BEG model there are two order 
parameters, the complete phase diagram is 3-dimensional, but we restrict ourselves to the 
2-dimensional sections just mentioned. In figure [I] we show the concentration/temperature 
(x-T) section and in figure ^| the chemical potential/temperature (A-T) section of the phase 
diagram. Both diagrams are obtained for a small value of the quadrupolar coupling constant 
K = 0.16. 

Next, we consider a special limit of the BEG-SK model to make contact to the work of 
Sherrington and Kirkpatrick. We allow the concentration to vary freely and we take the limit 
A — > — oo. From this, we recover the Sherrington and Kirkpatrick Hamiltonian all spins 
tend to be \s\ = 1, i.e., x — > 0. The phase diagram obtained by Sherrington and Kirkpatrick 
represents a section of the general phase diagram, i.e., the x-Jq-T section for a concentration 
value x = 0. This diagram exhibits a paramagnetic phase at high temperatures and small 
ferromagnetic coupling, whereas it exhibits a ferromagnetic phase for larger ferromagnetic 
coupling and moderate temperatures. Finally at very low temperatures the spin glass phase 
appears. This behaviour is shown in figure 0. The figure does not take into account the 
correction to the transition line between ferromagnetic and spin glass phase due to stability 
requirements M. 



The section of the general phase diagram obtained by Ghatak and Sherrington [|T^] is the 
A — T-section for Jo = 0, cf. ||22|| . The transition lines are obtained analytically following the 
procedure of Lage and de Almeida (T7| by examining the stability of the replica symmetric 



solution; they are given by 



A 



IT 



i±\/i-¥l+ Tln 



2(1-1) 

1 



T > - 

x ^ 3 



T < i 



(9) 



These transition lines are shown in figure f|. The lower two lines belong respectively to the 
two choices of sign in the second equation. 

Finally we are ready to present the complete phase diagram for BEG-SK model. Figure 
[5] shows the x-T section and figure ^| the A-T section of the phase diagram for one choice of 
the coupling parameters: Jo — 2.0, J = 1.0 and K = 0.16. In order to relate this to previous 
work and to present the complete phase diagram, the x-T and the A-T sections have been 
calculated for different values of the ferromagnetic coupling (scaled by the variance of the 
distribution of the couplings) and the results are shown in figures |7| and |8|, respectively. 
Figure [7| shows in the x-Jq/ J-T plane our numerical results (represented by diamonds) and 
the results of the SK model (represented by crosses) at a value x = of the concentration. 



The line of tricritical temperatures which will be given in formula (^) is shown with a 
broken line and the surface gives the second order phase transition. In figure |]the A-J /J- 
T section is given and the square symbols indicate the numerical results, which cover the 
first as well as the second order phase transition. The results from Ghatak and Sherrington 
|T3| are reproduced for the value of the ferromagnetic coupling J = 0, i.e., the face of the 
cube. For details on the stability we refer the reader to section [VI] and to figure [IT]. The first 
order transitions terminate at the line of tricritical potential given by formula (0), which 
is the thick line from the bottom left to the top right. The crosses represent the first order 



transition lines. The second order transition above T tri , given by formula (|2T| ) and obtained 
by the stability analysis presented in Section 6 are represented by squares. The stability 
lines below T tri are not given. 

We now discuss the numerical results in the light of recent work on disordered BEG-type 
models which are not frustrated. After we will address the problem of the stability of the 
replica symmetric solution. 



S 



V. BEHAVIOR OF THE FIRST ORDER TRANSITION 



Lets look at the influence of disorder on the first order transition, which can — following 
Blume, Emery and Griffiths — be interpreted as a phase separation. This first order transi- 
tion represents an example of a transition which breaks a global symmetry. The Hamiltonian 
is invariant under inversion of all spins; in the ordered phase this symmetry is spontaneously 
broken, whereas in the paramagnetic phase the symmetry is preserved. 

We first consider the non-disordered BEG model in the x — T representation, where 
the appearance of the phase separation is most visible. In figure [I] the phase diagram is 
shown for a weak quadrupolar interaction. In figure ^ the same is seen for a rather strong 
quadrupolar coupling K = 2.88. The point is that the quadrupolar coupling enhances the 
phase separation. In the limit of K » 1, the model tends to behave as the Griffiths model 
fl4| . We thus wish to emphasise the importance of the presence of a quadrupolar interaction. 
In the original model of Blume, Emery, and Griffiths, even the K = case exhibits a phase 
separation. It is due to the bulk interaction of the s = species (He 3 ) and the s = ±1 
species (He 4 ), the latter allowing an additional degree of freedom (ferromagnetic phase). At 
low temperature the ferromagnetic phase is the favoured phase for the s = ±1 species and 
phase separation occurs in order to permit ordering in a He 4 rich phase. Letting K ^ 
allows for an inter-isotopic interaction; in fact K represents a combination of inter-isotopic 
interactions, K = K33 + K44 — 2K34, which are assumed to be positive. It is remarkable that 
for a range of values of K the phase diagram exhibits a triple point together with the critical 
and tricritical point. This appears in the original BEG model and in the present model, as 



can be seen from figure |10[ The tricritical point is given in the diagrams by A, the triple 
point by B, and the critical point by C; these special points are also marked in figure 

Introducing disorder in the BEG model may affect the phase separation. Figure |5| shows 
the phase diagram obtained by a numerical extremisation of the free energy. The phase 
separation persists independently of the presence of disorder. The average over the Gaussian 
disorder has introduced an effective quadrupolar coupling k, which — as we have seen before 
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in the non-disordered case — enhances the phase separation. Furthermore this effective 
interaction is also temperature dependent: 

K = K+~@J 2 . (10) 

This promotes the phase separation at low temperatures as can be seen by comparing 
with figure |^ (for the low temperature regime) and figure p] (for the high temperature regime) 
in the non-disordered case. Even if there is an effect due to the presence of disorder, which 
might suppress or change the first order transition, this effect is apparently compensated by 
the additional term contributing to the effective quadrupolar coupling. 

This result complements those of Berker et al. ||, |§ and [|ITJ. They conjecture, in 
dimension d > 3, that the disorder lowers the tri critical temperature; then that part of the 
transition line enclosed between the former and the actual tricritical temperature should 
become second order. Also, they claim that all of the first order line is replaced by a second 
order one if the disorder is sufficiently strong. A few comments are in order: 

(i) The analysis of Berker et al. relies on the real space renormalisation group approach for 
an initially positive distribution of the couplings. This distribution remains positive 
and consequently cannot take into account the effect of frustration. Nevertheless, 
we do not think that the frustration is responsible for the persistence of the phase 
separation for the following reason: even when J is small, so that there is disorder 
but essentially no frustration, the phase separation persists (see the phase diagrams). 
Monte Carlo studies by Diep et al. || indicate this persistence and hence ensure that 
the phase separation is not an artefact of the mean field approach. 

(m) The mean field approximation is equivalent to a model in infinite dimensions, perhaps 
rendering the comparison of our results with the predictions of Berker et al. invalid. 
Nevertheless, it should give an indication of whether the first order transition per- 
sists or not. Furthermore one could argue that in infinite dimensions, an infinitely 
strong disorder is required to suppress the first order transition. However figure HI in- 



dicates that this model reduces to that discussed by Ghatak and Sherrington |13[ with 
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a vanishing ferromagnetic coupling, showing again the persistence of the first order 
transition. 

(Hi) The assumption of replica symmetry is not justified. However, based on the stabil- 
ity analysis, we have found that the first order transition persists; it seems rather 
improbable that it be suppressed by doing the full replica symmetry breaking scheme. 

To summarise, Berker et al. claim that their results are generic, but our mean field 
treatment shows that the phase separation is not affected by the introduction of disorder. In 
particular, the tricritical temperature is not modified: T / J = 1/3, if J > J , using the same 
scaling as Berker et al. Note that their analysis does not reveal the effective quadrupolar 
interaction, which plays a subtle role. On the one hand, following Berker et al. and their 
general arguments 011, the disorder should change the first order transition to a second order 
one; on the other hand, the effective quadrupolar interaction enhances the phase separation, 
so stabilises the first order transition. As we have seen above, the effective quadrupolar 
interaction more than compensates the first effect so that the first order transition persists. 
Furthermore, since for a small amount of disorder the frustration effects are negligible, the 
models considered are comparable and our results provide a counterexample to their claim. 

VI. TRANSITION LINES FROM A STABILITY ANALYSIS 

The assumption of replica symmetry allows for a simple answer for the BEG-SK model. 
However, as for the SK model, it leads to an unstable solution and the breaking of the 
permutation symmetry of the replicas is required. We expect the Parisi breaking scheme 
to apply here and to remove the instabilities in our model. Since under replica symmetry 
breaking the findings of the analysis assuming replica symmetry are in general confirmed, 
i.e., a spin glass phase persists, we won't study finite nor infinite step step replica symmetry 
breaking. In this section we restrict ourselves to the replica symmetric solution and we 
derive the lines of instability to locate the phase transitions lines given in section IV — even 
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though the location of these lines may differ in a more advanced treatment. The stability 
of the replica symmetric solution of the SK model in its integer spin generalisation was 
examined by Lage and de Almeida [|T7|] and in greater detail by Mottishaw and Sherrington 



p9| . These works limit their analysis to the disordered integer S spin glass model without 
a non-zero ferromagnetic mean and without quadrupolar coupling. 

The present stability analysis follows closely the moethods of de Almeida and Thouless 
and subsequent works. To examine the stability of the solution it has to be made sure 
that the solution extremises the free energy. The quadratic form T describing the deviation 
of the solution from its stationary value should be positive definite. 

(3f = (3f(m,z,q)-^T + 0(5 3 ) and 

a/3 a/3 (a/3)(rfS) 

a/3 (a/3)-y (a/3)7 

The matrix G associated with this form is the Hessian TC. The eigenvalues of the Hessian 
should be non-negative to make sure that the solution is stable. Due to the same symmetry 
arguments used by de Almeida and Thouless J/J], we discover three families of eigenvalues 
which by convention are called: longitudinal, longitudinal anomalous and replicon. But in 
the present case the stability analysis requires the diagonalisation of a Hessian 7i which 
is built out of three blocks. As a consequence there are now three distinct eigenvalues in 
each family. The longitudinal eigenvectors are the fully symmetric ones with respect to 
permutations of the replica indices and are of the form: 

= a for fj, — 1 . . .n, — b for v — 1 . . .n, e^g = c for a, j3 — 1 . . . — (n — 1). (12) 

They give rise to eigenvalues which are the solution to the cubic characteristic equation. 
The full formula for the eigenvalues is too cumbersome to be displayed here, they will be 
discussed in the paramagnetic phase below (the full expressions are given in appendix 
The longitudinal anomalous eigenvectors are generated by one distinct replica index 9 and 
look like: 
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= a for /x = 1 . . . n and /x 7^ = g for /i = 9, 

e L v A = b for v = 1 . . . n and /x ^ 9, e L v A = c for v = 6, (13) 
e a@ = d for a, /3 = 1 . . . — (n — 1) and a, (3 ^ 9, e^p = e for a, /3 = 0. 

Again the eigenvalues are given in appendix The replicon eigenvectors are generated 
by two distinct replica indices u and 9 and they are of the following form: 
= a for /x = 1 . . . n and /x 7^ 6 1 , u;, = g for /x = 6 1 , a;, 
= b for z/ = 1 . . . n and /x 7^ 9,uj,e^ = c for u = 9,lu, 

(14) 

e f/3 = for a, /5 = 1 . . . — (n — 1) and a, f3 ^ 9 and a, (3 ^ u, 
e?g = e for a, /3 = or a, (3 = to, efg = f for (a(3) = (9tu). 

The eigenvalue is given by the solution of the characteristic equation, but as this equation 
is explicit, so the eigenvalue is obtained immediately: 

X R = P-2Q + R, (15) 

where P, Q and R are appropriate correlation functions given in the appendix. Adding the 
number of eigenvalues obtained in the three symmetry families gives the expected number 
3 + (3rx — 3) + |(rx — 3) = |(n + 3), i.e., the dimension of TC. The so called longitudinal eigen- 
values Xl and Xla coincide in the n — > 0-limit as in the work of de Almeida and Thouless 
and Lage and de Almeida [[L7j]. The replica symmetric fluctuations of the one- indexed 
quantities (z a and m a ) and the two- indexed quantity q a p are described respectively by the 
longitudinal eigenvalues, Xla and Xl. The third eigenvalue Xr is distinct from the other 
two eigenvalues and is related to the fluctuations which break the replica symmetry of the 
two indexed quantity g Q( g. In order to investigate the stability of the replica symmetric 
solution, for example in the ferromagnetic phase, the eigenvalues Xl, Xla and Xr have to 
be calculated in this phase. In the appendix the general expressions for the different eigen- 
values and the correlation functions appearing therein are given. They must be evaluated 
under the appropriate approximation, i.e., in the present case, under the replica symmetric 
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approximation. In order for a phase to be stable the eigenvalues in this phase have to be 
non-negative. We now investigate the stability of the paramagnetic phase. 



The paramagnetic phase 

Consider the paramagnetic phase defined by the vanishing of all magnetic order pa- 
rameters: q = 0, m = 0. This considerably simplifies the expressions for the eigenvalues. 
All off-diagonal elements of the Hessian or combinations appearing in the characteristic 
equations vanish. In order to simplify the notation further, we consider directly the combi- 
nations appearing in the definitions of the different eigenvalues. The longitudinal eigenvalues 
are given by the solution of: 

A| + a\ 2 L + b\ L + c=(X- Ai)(A - A 2 )(A - A 3 ) = 0. (16) 

The coefficients for a cubic equation are given implicitly by the solutions to that equation: 

- a = [U - V] + [A - B] + [P - 4 Q + 3 R] = Ax + A 2 + A 3 

b=[U- V\[P -4Q + 3R] + [A- B\[P - 4Q + 3 R] + [U - V][A - B] = A : A 2 + A 2 A 3 + AxA 3 
c= -[U-V][A-B][P-AQ + 3R] = -AiA 2 A 3 , 

which may be solved immediately to give 

Xi = [A-B] 

\ L 2 =[U- V] (18) 
A3 = [P-4Q + 3R] = P. 

The longitudinal anomalous eigenvalues are given by the same expressions. The replicon 
eigenvalue is given by: 

\ R = P-2Q + R = P, (19) 

and coincides in the paramagnetic phase with the third longitudinal eigenvalue. In order for 
the paramagnetic phase to be stable, the eigenvalues in this phase must be non-negative. 
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The border line of stability is given by \l,la,r — 0, which defines seven lines; but due to 
the collapsing eigenvalues of different families there are only three different stability lines. 
In the replica symmetric approximation, these read: 



A-B = (3J [1 -/3Jqz] = 

U — V — 0k[1 - [3k(z - z 2 )} = 



(20) 



P-AQ + 3R = P-2Q + R = 2 J 2 [l - (3 2 J 2 z 21 



0. 



The stability limits are given by the following explicit formulae, first in the plane of 
concentration x versus temperature T: 

T(x) 



Jo 
T(x) 



J 



1-x 

1 — X 
x(l — x). 



(21) 



The paramagnetic phase becomes unstable whenever the temperature is lowered below 
the greater of the two bilinear couplings J and Jo- This behaviour is identical to that 
observed by Sherrington and Kirkpatrick. The conditions of stability can be rewritten in 
terms of the chemical potential A and temperature T. The stability lines are obtained by 
using the fixed point equation (H) determinating the concentration, and represent a line in 
the A-T-section of the phase diagram. 



A R s(T) 



(K + lf3J 2 )^ + T\n [2 (^-l) 



{K + \p.J 2 )\(l±J 



1 - 



2f3K+/3 2 J 2 



) +Tln 
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T7I, 



T > T tri 
T < Ttrf 

(22) 



For one choice of the parameter the complete set of stability lines has been depicted in 
figure Furthermore the numerical results covering the first order transition points as 
well as those of second order are represented by diamond symbols. The following discussion 
refers to this figure. The two choices of sign in the second equation separate different regions 
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of the free energy's behaviour. The branch (/i) belonging to the positive sign separates 
the region of free energy with a unique minimum (PM-U) from the region where there is 
more than one local minima (PM-M). It is in this latter region where first order transitions 
occur and the phase separation appears. The branch (g) belonging to the negative sign 
represents regions where the free energy has three local minima, separating the region where 
the minimal solution is the one with vanishing order parameter from the region where the 
minimal solution is the one with a nonvanishing order parameter; this line (g) is indicating 
the onset of phase separation. The concentration x jumps when crossing this line; in fact 
this line represents in the A-T section of the phase diagram the whole phase separation or 
coexistence region. As can be seen by comparision with the numerical data, the location 
of the first order phase transition is predicted correctly only near the tricritical point. It is 
known that the stability analysis is not an adequate mean to determine first order phase 
transitions. As in previous works the A-line (A) meets the first order transition line at a 



tricritical point: a second order line changes to a first order line, or, following Griffiths [15 



three critical lines meet (see also [0). For completeness the continuation of the A-line below 
T tri has been drawn too. 



Tricriticality 

The occurrence of a tricritical point in a system with one order parameter is signaled by 
the vanishing of the fourth derivative of the free energy with respect to the order parameter. 
In the present problem there are three order parameters and they are have replica indices. 
The criterion must be modified: the eigenvalues of the matrix of the fourth derivatives have 
to vanish. But due to the diagonal structure of the Hessian TC in the paramagnetic phase 
this amounts to computing the fourth derivative with respect to the magnetisation and to 
search for its zero. This is in fact the classical argument of the Landau theory. The tricritical 
temperature in the non-disordered case is given for a vanishing quadrupolar coupling by 

T tri = i max {J , J} . 
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For a non-vanishing quadrupolar coupling the tricritical temperature is given by: 

The tricritical point present in the non-disordered model (BEG) persists and in terms of 
the temperature and concentration is given in the crudest approximation by: 

£ t ri = \ and T tri = ^ + \ max { J , J} . (23) 
6 ZK + 6 

The numerical data show that this is a good approximation. It is worth mentioning 
that the tricritical concentration has kept its value from the non-disordered BEG model. 
In contrast the tricritical chemical potential has been modified. Resolving the tricritical 
condition provides a formula in A-Jq/J-T space for the tricritical line. In the A-T plane, it 
is given by: 

\ f + 57 + \ J o ln[2] and T tri = f for J > J 
A tri ( J , J) = <( 3 . (24) 

f + \J+ | J In [2] and T tri = f for J > J 

This line, as a function of Jo, is constant up to Jo — J, and then tends for stronger 
ferromagnetic coupling to larger chemical potential. The line of stability for the second 
order transition and the aforementioned line of tricritical A are displayed together with 
the numerical data in figure |8|. The value for the quadrupolar coupling for this figure is 
K = 0.16. 



VII. DISCUSSION AND CONCLUSIONS 

We have performed a replica study of the disordered BEG model and have extended 
previous work to present a picture of the complete phase diagram. The BEG-SK model 
shows a second order transition, the so-called A-line, separating an ordered phase from 
a paramagnetic one. The ordered phase may be a ferromagnetic or a spin glass phase, 
depending on the strength of the ferromagnetic coupling Jo. Furthermore, a transition of 
first order, present in the non-disordered BEG model, persists, and may be interpreted as 
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a phase separation. This extends recent work of Berker et al. on the influence of disorder 
on first order transitions. Our mean field study is completed by a stability analysis of 
the replica symmetric approximation. The complete set of eigenvectors and eigenvalues 
has been found and analysed in the paramagnetic phase. The replica symmetric solution 
suffers from instabilities, but the assumption of replica symmetry and the associated results 
do not exhibit more severe inconsistencies than in the SK model. For the SK model, the 
instabilities of the replica symmetric solution have been remedied by the infinite step replica 
symmetry breaking procedure. We expect the same approach to work in this case, but 
following Mottishaw and Sherrington near the tricritical point, the Parisi Ansatz has to 
be extended to higher order in q a p. 
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APPENDIX A: EIGENVECTORS AND EIGENVALUES 



We follow Lage and de Almeida |Ty|] in their stability analysis of the replica symmetric 
solution. The Hessian of the free energy is given schematically by: 
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n 



d 2 f 



d[z a m a q a p}d[z a m a q Q p} 



U{z,z) V 
U 



W X 

w 

A(m,m) B 
A 



V 



Y Z 

Y 
C D 

C 

P(q,q) Q R 
P Q 
P 



The second equality defines the quantities U, V, . . . with respect to their position in the 
Hessian. These quantities are the respective Gff , etc ... as given in the main part. As will 
be seen later they can be expressed in terms of the different multispin correlation functions 
up to degree 4. The dimension of the Hessian is n + n + — 1) = |(n + 3) and equals 
the number of eigenvalues and eigenvectors to be found. In order to diagonalise the Hessian 
we construct the eigenvectors using the symmetry arguments exposed by de Almeida and 
Thouless . We do not construct a orthonormal set of eigenvectors, because the additional 
constraints destroy the symmetry with respect to replica permutations. Rather we content 
ourselves with three families of eigenvectors, each orthogonal to another, but not orthogonal 
within the families. 

We start with the eigenvector totally symmetric under permutations of the replicas: 

a : /i = 1 . . . n 
£l = { b : v = 1 . . . n 

c : a,P= l...f(n- 1) 

According to conventional notation this vector is referred to as belonging to the longi- 
tudinal subspace; hence the subscript L. This is a subspace of dimension d = 3, which is 
easily verified by constructing a orthogonal set of eigenvectors conserving the form prescribed 
above. The eigenvector equation TCeL = Xl^l reads 

a[U + (n- l)V - \ L ] + b[W + (n - 1)X] + c[Y(n - 1) + Z[-(n - 1) - (n - 1)]] = 



19 



a[W + (n- 1)X] + b[A + (n - l)B - \ L ] + c[C(n - 1) + D[-(n - 1) - (ra - 1)]] = 

72 

a[2F + (n - 2)Z] + b[2C + (n - 2)D] + c[P + 2(n - 2)Q + i2[-(n - 1) - 2(n - 1) - 1] - A L ] = 

The longitudinal eigenvector equation gives rise to a cubic characteristic equation for the 
eigenvalue A^: 

(A L ) 3 + a(A L ) 2 + 6A L + c = 0. 

We give its coefficients in the n = limit: 

a = —[U -V] - [A - B] - [P-4Q + 3R] = -\[ - \% - A3 
6 = -[1^ - X] 2 + [U - V][P - 4Q + 3 R] + 2[Y - Z} 2 

+[A - B\[P - 4Q + ZR] + [U - V\[A -B} + 2[C - Df = Af A^ + A^Af + AfAf 
c = -[17- V][A - B][P - 4Q + 3R] - 2[U - V][C - D} 2 

+ [W - X] 2 [P -4Q + 3R] + 4[W - X][Y - Z][C - D] - 2[Y - Z] 2 [A - B] = -A^A^Af. 



Define the following quantities: 



° 2 , 

7 = ~¥ 

/a\ 3 1 , 



37 V2 

where only the real cubic roots are used. This allows us to write the three solutions as: 
A L i = <j + + <j~ 

\ L2 = - l -(a + + a') + ^V3(a + - a~) (Al) 
A« = -±(<7 + + 0-^(<7 + -0- 
The value of V tells us whether the solutions are degenerate or not. 
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< : distinct real solutions 
= : degenerate real solutions 
> : complex solutions 



As has been recognized by Lage and de Almeida |17|], occasionally T is positive and so 
the eigenvalues become complex. This has been verified numerically by da Costa et al. 0. 

The next eigenvectors, which will be called longitudinal anomalous, are constructed by 
breaking the symmetry of the longitudinal vector with respect to one replica, given by the 
distinct index 9. 

fi = 1 . . .n, fjb ^ 9 
li = 9 

v = 1 . . . n, /i ^ 9 

v = 9 

a,P = l...*(n-l),a,P^e 
a, = 9 

Consider the orthogonality condition for the first eigenvector, if a = g or b = c. This 
condition results in trivial eigenvectors. In order to obtain non-trivial eigenvectors, the 
symmetry in both n-blocks has to be broken, giving kn eigenvectors of the second family, 
where k is the number of different choices of the parameters conserving the prescribed form. 
Orthogonality of the second family to the first family of eigenvectors requires: 



9 
b 

c 
d 
e 



g = (1 — n)a c = (1 — n)b 



(A2) 



Therefore k = 3 and this choice gives rise to 3n eigenvectors including the previous one. 
Writing down the characteristic equation for the longitudinal anomalous eigenvalues Xla 
gives: 

ii 

a[U — V — \ LA ] + b[W -X} + d{- -1)[Y-Z]=0 
a[W -X} + b[A-B- X LA ] + d(- — 1)[C — D] = 

T) 

2a[Y -Z} + 2b[C — D] + d{- - l)[P + (n - 4)Q + R(3 - n) - X LA ] = 
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The longitudinal anomalous eigenvector equation also gives rise to a cubic characteristic 
equation for the eigenvalue A^: 

(A LA ) 3 + a(X LA ) 2 + b\ LA + c = 0. 

We quote its coefficients in the n = limit: 

a =-U+V-A+B-P+4Q-3R 

b = (W - X) 2 - (U - V) (P - 4 Q + 3 R) + (2 Y - 2 Z) (-Y + Z) 

- (A - B) (P - 4 g + 3 R) - (U - V) (A - B) + (-C + D) (2 C - 2 D) 
c = (U - V) (A - B) (P - 4 Q + 3 R) - (U - V) (-C + D) (2 C - 2 D) 

- (W - X) 2 (P - 4 Q + 3 i2) + (W - X) (-Y + Z) (2 C - 2 D) 

+ (2Y -2Z)(W - X) (-C + D) - (2 F - 2 Z) (-F + Z)(A- B) 



Again the solutions may be written as in equation |A1|, but now with modified coefficients 



7,£, .... 

The eigenvalues will in general be different in each family, but since the characteristic 
equations for the longitudinal and the anomalous longitudinal become identical in n = 0- 
limit, so do the eigenvalues, and then the longitudinal and longitudinal anomalous families 
collapse. There are 3n eigenvectors and the three eigenvalues are each n-fold degenerate, 
including the three eigenvalues of the first eigenvector. 

Finally, there is the third family, called replicon. Breaking the symmetry with respect to 
permutations of pairs by distinguishing two indices 6 and lu, we obtain the following form 
for the eigenvectors: 
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a 


fi= 1 




n, fj, 7^ 9, uj 


g 


fi = 9, 


uj 




b 


v = 1 




n, [x 7^ 9,uj 




v — 9 

IS — u . 






d 


a, (3 — 


1 


. . . |(n — 1), a, ft ^ 9 and a, j3 ^ uj 


e 


a, /3 = 


9 


or a, (3 = uj 


f 


(aP)~- 


= ( 


9uj) 



Orthogonality of the third family to the previous two families of eigenvectors requires: 

g = a c = b = /+(n-4)e+(3-n)d. (A3) 

The choice of equal non- vanishing entries would give the first family, so the recommended 
choice is 

g = a = c = b = /=(2-n)e e=~(3-n)d (A4) 
The characteristic equation is the solution for the replicon eigenvalue itself: 

\r = P — 2Q + R. 

The replicon eigenvalue is independent of n and identical in form to the result obtained 
by de Almeida and Thouless 0. This eigenvalue is § (n — 3)-fold degenerate. 

1. The replica symmetric Ansatz 

In order to analyse the stability of the replica symmetric Ansatz, the quantities defined 
in the Hessian are to be evaluated assuming replica symmetry. This gives: 

A = (3J [1 - PU< s 2 a >-<s a > 2 )] = /3J [1 - pj (z - m 2 )] 

B = P 2 Jl[< s a > 2 -< s a s p >} = l3 2 Jl\m 2 - q] 

C = f3 2 J 2 /3J [< s a >< s a S/3 > - < s a s 2 p >} = (3 2 J 2 /3J [mq - r] 
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D = p 2 J 2 pJo[< s 7 >< s a Sf3 > — < SaSpS-f >} = p 2 J 2 pJo[mq — u] 

P = P 2 J 2 [1 - (3 2 J 2 (< sis 2 , >-< s a s p > 2 )] = P 2 J 2 [1 - P 2 J 2 (s - q 2 )] 

Q = ^j 4 [< s a Sp > 2 -< s 2 a sps, >} = p 4 J 4 [q 2 - v] 

R = f3 A J A [< s a sp > 2 - < SaSps^ss >} = (3 A J A [q 2 - w] (A5) 
U = Pk[\ - (3k(< s 2 q > - < s 2 a > 2 )} = (3k[1 - (3k(z - z 2 )} 

V = /?V[< S 2 a > 2 -< S 2 a sl >] = P 2 K 2 [Z 2 - S] 

W = PJ Pk[< s a > (< s 2 a > -1)] = pJ pK[m(z - 1)] 

X = (3Jq(3k[< s a >< s 2 a > — < s a s 2 p >} = j3Joj3K[mz — r] 

Y = (3 2 J 2 /3k < s a s p >[<s 2 a > -1] = p 2 J 2 {3Kq[z - 1] 

Z = /3 2 J 2 /3k[< s 2 a >< s a Sf3 > - < s 2 a s/3S y >} = (3 2 J 2 (3K,[zq - v]. 

The second equality in each line results from the assumption of replica symmetry. Addi- 
tional simplifications are due to the fact that for the S = 1 spin model there are additional 
relations, e.g., S 2 = S 4 , etc. Using the previously defined functions 4>k(y), see equation (0), 
the replica symmetric correlation functions introduced above may be explicitly written as: 



m = 


<s a >= J 


v vMv) 


t = 


< slsp > 


z = 


<sl >=/ 


v vMv) 


u = 




q = 


< S a S{3 > = 


J v yUv? 


V = 


< S^SpS^ 


r = 


< S a Sp > = 


J v yMv)My) 


w = 




s = 


< 44 >= 


J v vMy) 2 







VyMy? 

VyMvYMv) (A6) 

VyUyf- 



2. Paramagnetic phase 

As is seen from the definitions of the multispin correlation functions (see equations [AS 



and |A6|) some of these vanish identically in the paramagnetic phase. Furthermore some of 
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the combinations appearing in the calculation of the eigenvalues, e.g. [W — X] = 0, vanish 
too. This simplifies significantly the stability analysis in the paramagnetic phase, because 
the Hessian Ti and the matrix of the fourth derivatives become diagonal. 
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FIGURES 




x 

FIG. 1. The phase diagram of the (non-disordered) BEG model in the concentra- 
tion/temperature-plane, i.e., x — T-plane, for the bilinear coupling Jo = 1.0 and for weak quadrupo- 
lar coupling K = 0.16. The thin line represents a transition of second order, whereas the diamonds 
indicate a transition of first order. FM/SG stands for ferromagnetic/spin glass phase,, PM is the 
paramagnetic phase; the point A indicates the tricritical point. 
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FIG. 2. The phase diagram of the (non-disordered) BEG model in the chemical poten- 
tial/temperature-plane, i.e., A — T-plane, for the bilinear coupling Jo = 1.0 and for weak quadrupo- 
lar coupling K = 0.16. The thin line represents a transition of second order, whereas the diamonds 
indicate a transition of first order. 



T/J 




0.2 0.4 0.6 0.8 
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Jo/J 

FIG. 3. The phase diagram of the SK model in the ferromagnetic coupling/temperature-plane, 
i.e., Jq/J — T-plane. The lines shown are transitions of second order. 




FIG. 4. The phase boundary of the GS-model (K = 0.0 and Jo = 0.0) in the chemical poten- 
tial/temperature-plane, i.e., A — T-plane. 
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FIG. 5. The phase diagram of the disordered BEG-SK model in the x — T-plane for the bilinear 
couplings Jo = 2.0, J = 1.0 and for weak quadrupolar coupling, K = 0.16. The thin line and the 
dense lying diamond symbols represent the transition of second order, whereas the sparse lying 
diamonds indicate the first order transition. The two dotted lines are the free energy minima for 
two different chemical potentials (A = 0.6 and 0.65) and varying the temperature (from T = 0.0 
to T = 1.2). 
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FIG. 6. The phase diagram of the disordered BEG-SK model in the A-T-section for the bilinear 
couplings Jo = 2.0, J = 1.0 and for weak quadrupolar coupling K = 0.16. The thin line and the 
dense lying diamond symbols represent the second order transition, whereas the sparse diamonds 
indicate the transition of first order. 




FIG. 7. The complete phase diagram of the disordered BEG-SK model, with the lines 
of stability and the line of tricritical temperature in the concentration/ferromagnetic cou- 
pling/temperature-representation, i.e., x— Jo/ J— T-plane, for weak quadrupolar coupling K = 0.16. 
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FIG. 8. The complete phase diagram of the disordered BEG-SK model, with the lines of 
stability and the line of tricritical chemical potential in the chemical potential/ferromagnetic 
coupling/temperature-representation, i.e., A — Jo /J — T-plane, for weak quadrupolar coupling 
K = 0.16. 
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FIG. 9. The phase diagram of the non-disordered BEG model in the x — T-plane for the bilinear 
coupling J = 1.0 and for a rather strong quadrupolar coupling, K = 2.88. The point A indicates 
the tricritical point, the point B represents the triple point and point C is the critical end point. 
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FIG. 10. The phase diagram of the non-disordered BEG model in the A — T-section for the 
bilinear coupling J = 1.0 and for a rather strong quadrupolar coupling K = 2.88. 
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FIG. 11. The phase diagram of the disordered BEG-SK model, with the lines of stability 
in the chemical potential/temperature-representation, i.e., A — T-plane, for weak quadrupolar 
coupling, K = 0.16, and for the ferromagnetic coupling Jo/ J = 2.0. The symbols have the 
following meanings: PM-M is the paramagnetic phase with many free energy minima, PM-U is the 
paramagnetic phase with an unique free energy minimum, FM/SG the ferromagnetic or eventually 
spin glas phase, \i gives the positive sign branch and g the negative sign branch of equation (|22|); 
A is the second order transition line. 
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